Mon. Not. R. Astron. Soc. 000, [|-|lo| (0000) Printed 5 February 2008 (MN I*TeX style file vl.4) 



The initial conditions of isolated star formation: IV — C O 
observations and modelling of the pre-stellar core L1689B 



N. E. Jessop ^ and D. Ward-Thompson^ 

^ Joint Astronomy Center, 660, N. A 'Ohoku Place, Hilo, Hawaii, USA 

^ Department of Physics and Astronomy, Cardiff University, P.O.Box 913, Cardiff, UK 



Accepted 2000 October 25. Received 2000 October 24; in original form 2000 April 19. 



o 
o 
o 

(N 
o 

CD 

Q 
in 



> 
m 
a^ 
o 

(N 

o 
o 

o 



.;-( 1 INTRODUCTION 



ABSTRACT 

We present C^^O observations of the pre-stellar core L1689B, in the (J=3— >2) and 
(J=2— >!) rotational transitions, taken at the James Clerk Maxwell Telescope in 
Hawaii. We use a A-iteration radiative transfer code to model the data. We adopt 
a similar form of radial density profile to that which we have found in all pre-stellar 
cores, with a 'flat' inner profile, steepening towards the edge, but we make the gradient 
of the 'flat' region a free parameter. We find that the core is close to virial equilibrium, 
but there is tentative evidence for core contraction. We allow the temperature to vary 
with a power-law form and find we can consistently fit all of the CO data with an in- 
verse temperature gradient that is warmer at the edge than the centre. However, when 
we combine the CO data with the previously published millimetre data we fail to find 
a simultaneous fit to both data-sets without additionally allowing the CO abundance 
to decrease towards the centre. This effect has been observed qualitatively many times 
before, as the CO freezes out onto the dust grains at high densities, but we quantify 
the effect. Hence we show that the combination of mm/submm continuum and spec- 
tral line data is a very powerful method of constraining the physical parameters of 
cores on the verge of forming stars. 

Key words: stars: formation - ISM: globules. 
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Star formation occurs in dense molecular cloud cores, and 
many surveys of such regions have previously been carried 
out, including the pioneering work of Myers and co-workers 
(e.g. Benson & Myers 1989 and references therein). They 
separated these cloud cores into those that had already 
formed stars and thus contain embedded Young Stellar Ob- 
jects (YSOs), and those that had not - the so-called 'starless 
cores' (Beichman et al 1986). The starless cores are prime 
candidates to study observationally the sites of potential fu- 
ture star formation, as they are believed to represent the 
initial conditions for protostellar collapse. We have been ob- 
serving starless cores for a number of years to try to con- 
strain theoretical models of protostellar collapse. 

Ward-Thompson et al. (1994 - hereafter Paper I) 
showed that many starless cores contain dense central con- 
densations which they named 'pre-protostellar cores' (or 
more recently 'pre-stellar cores' for brevity). Detailed ob- 
servational studies of pre-stellar cores offer the opportunity 
to ascertain the density and temperature distribution within 
the core, as well as the kinematics and details of the chem- 
istry, including dust-molecule interactions. AH of these fac- 



tors are thought to play important roles in governing the 
way in which a body of gas collapses to form a protostar. 

In Paper I we found that the variation of density (p) 
with radius (r) in pre-stellar cores is very different from the 
singular isothermal sphere (p oc everywhere) originally 
suggested by Shu (1977) as the initial conditions for star for- 
mation. Instead the cores appear to have a much flatter den- 
sity profile in the inner region (p cx r~^), steepening towards 
their edges (p cx f~^). This was subsequently confirmed for 
the pre-stellar core L1689B by Andre, Ward- Thompson & 
Motte (1996 - hereafter Paper II), and for other cores by 
Ward- Thompson, Motte & Andre (1999 - hereafter Paper 
III). Most recently, an ISOCAM study by Bacmann et al. 
(2000) has shown that some cores have very steep edges in- 
deed (p cx r~^). In a study of one prestellar core (L1544), 
Tafalla et al. (1998) discovered significant large-scale mo- 
tions of gas, possibly indicating contraction of the core. 

In this paper we present a C^^O study of one pre-stellar 
core, L1689B, which has been previously studied at 1.3mm 
(Paper II), and by using a detailed model of the radia- 
tive transfer, we investigate the physical parameters of this 
core. The remainder of the paper is laid out as follows: Sec- 
tion 2 presents the C^*0 observations of L1689B. Section 3 
presents a radiative transfer model of L1689B. In Section 4 
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R.A. offset(arcsec) from 16 31 46.98 

Figure 1. L1689B as seen in the C^**0 (J=2^1) rotational transition. Tmb/K is plotted against v/kms~^ at each position. Axes are 
labelled in 1950 co-ordinates. See text for discussion. 



the variation of density and temperature of the gas is ex- 
plored, and we introduce the 1.3mm data to show how that 
leads to further constraints, including the variation of the 
abundance of CO relative to the dust. Section 5 summarises 
the main conclusions. 



2 C^^^O DATA 
2.1 Observations 

The observations were carried out at the James Clerk 
Maxwell Telescope (JCMT)0, located on Mauna Kea, 
Hawaii, on 1995 July 18th, 22nd and 25th 17:30-01:30 HST 
(UT = 03:30-11:30), on 1996 August 31st and September 
1st 17:30-01:30 HST (UT = 03:30-11:30). 

The C'^^O (J=3^2) transition, with a rest frequency 
of 329.33 GHz, and the C^**0 (J=2^1) transition, with a 
rest frequency of 219.56 GHz, were observed using the com- 
mon user heterodyne receivers RxB3i (Cuimingham et al. 



1992) and RxA2 (Davies et al. 1992). The JCMT half-power 
beam-width (HPBW) is 19 arcsec at 220 GHz and 14 arc- 
sec at 329 GHz. Double-sideband system temperatures were 
2000-12000 K for receiver B3i and 480K for receiver A2 ob- 
servations. The backend used was a digital auto-correlation 
spectrometer (DAS), with a resolution of 378kHz and 95 
kHz per channel, for the J=(3— >2) and J=(2^1) transitions 
respectively, corresponding to 0.34 kms~^ at J=(3^2) and 
0.13 kms~^ at J=(2— *1). Data reduction was carried out us- 
ing the SPECX package (Padman 1990). The weather during 
the observations was in general good for millimetre (RxA2) 
observations and short periods proved adequate for obser- 
vations of the submillimetre J=(3^2) transition. 

Regular pointing checks on bright sources, observed 
with the continuum backend system, were carried out 
throughout the observing run to check the performance of 
the telescope. We found that at worst the pointing was accu- 
rate to within ~3 arcsec. Once the receivers had been tuned 
to the correct frequency, a bright standard source was ob- 
served to check the chopper wheel calibration. 



* JCMT is operated by the Joint Astronomy Center, Hawaii, 
on behalf of the UK PPARC, the Netherlands NWO, and the 
Canadian NRC. 
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2.2 Data Reduction 

The data obtained were in the form of a series of spectra 
sampled at different positions that could be built into a data 
cube. The basic data format is a spectrum whose intensity 
is given in terms of the antenna temperature T^. The spec- 
tra were calibrated regularly using the usual chopper wheel 
method (Kutner & Ulich 1981) in addition to the standard 
observations described above. 

Since the beam of the telescope is not expected to be 
perfectly coupled to the source, a correction factor 77 has to 
be applied to T* in order to calculate the radiation tempera- 
ture of the source Tmb- This coupling constant is dependent 
on the spatial extent of the source. Taking into account the 
observations of standard sources discussed above and using 
the suggested coupling factor of Matthews (1999), t^b, for 
a source filling the beam, we used the following conversions 
from antenna temperature to main beam temperature for 
receiver A2: T^b = T*/0.65(±5%); and for receiver B3i: 
Tmb = T*/0.4(±10%). These conversion factors are slightly 
lower than the canonical values, which may have been due 
to minor receiver problems at the time of the observations. 

2.3 L1689B C^**0 Data 

L1689B data cubes in the C^*0 (J=2^1) and C^^O 
(J=3^2) transitions are presented in Figures 1 and 2. The 
source is clearly detected in both transitions at virtually 
all positions in the maps. The brightness distribution of the 
core appears relatively uniform across the area mapped, and 
it is not strongly peaked. The noise level per channel in 
the spectra is 0.6 K for C^^O (J=2^1) and 0.4 for C^**0 



(J=3^2). For most positions Tn,b(peak) for C^*0(J=2^1) 
lies between 5 and 6 K. The C^*0(J=3^2) T^b intensity 
at most positions lies in the range 3.8 to 4.8 K. The most 
clearly detected variation in intensity seems to be a drop off 
towards the North-East of both maps, most noticeably in 
C'«0(J=3^2). 

The full-width at half maximum (FWHM) width of the 
lines is 0.6-0.8km/s. The C^**0(J=3^2) map had too little 
signal to noise to reveal any kinematic information, but a 
channel map of the C^*0(J=2— >1) map was produced and 
is shown in Figure 3. The map shows no evidence of rotation 
which would be revealed by the red and blue shifted maps 
peaking at opposite positions either side of the centre of the 
core. However some evidence exists for the points 20 arcsec- 
onds east and west of centre being blue shifted with respect 
to the central position. Detailed inspection of the individual 
lines shows this to be true, as they are shifted by approx- 
imately 0.2-0.3 km/s. However, there is no clear evidence 
for line asymmetry. We hypothesise that this may be due to 
slow contraction of the core, but this would require further 
observations to confirm. Using the equation for the virial 
mass, A'fvir, of a spherical cloud (e.g. MacLaren, Richardson 
and Wolfendale 1988): 



where v is the line width in km/s, R is the radius in pc and k 
is a constant (between 126 and 210, depending on the exact 
form of the density distribution), we derive a virial mass 
for the core of <6-8Mq. In Paper II we estimated the mass 
within 120 arcsec to be >2.1M0, assuming a temperature 
of 18K. This lower limit may be an underestimate since in 
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Figure 3. A channel map of the L1689B in the C^^O (J=2— >1) transition. Each map in the sequence shows the integrated intensity of 
a 0.2km/s range in velocity. The line centre is at 3.7 km/s. The red-shifted maps peak on centre. However the blue shifted maps peak to 
the east and west, due to a slight shifting of the position of the line centre at these positions. No evidence for rotation exists. The scale 
bar shows the gray scale in units of K. 



Paper II we used a temperature of 18K, while more recent 
studies (e.g. Bacmann et al. 2000) found a temperature of 
~12K. This would give a mass of ~4.2Mq, closer to our 
virial estimate. Hence we see that the core is close to virial 
equilibrium. 

To quantify the spatial structure in the maps we chose 
to azimuthally average the spectra in the C^*0 (J=3— »2) 
and C^*0 (J=2— >1) maps, centering the averaging around 
the C^^O (J=2— >1) peak, which is coincident with the mm- 
continuum peak. In terms of angular displacement from this 
position, the maps contained in each transition: 2 spectra at 
7 arcsec; 4 at 16 arcsec; and 2 each at 22, 26 and 32 arcsec 
displacement from the centre. 

We can summarize the results of the azimuthal aver- 
aging by saying that; the C^*0(J=2— >1) Tjnb(peak) at the 
centre is 5.5K; the ratio of C^*0(J=2^1) at a distance of 
32 arcsec from the centre, to C^*0(J=2^1) at centre is 
IzbO.l (much less centrally peaked than the 1.3mm contin- 
uum emission) ; the ratio ofC^*0(J=3^2)toC^^O(J=2^1) 
is 0.78 ± 0.09 at centre; the line width is 0.7km/s ±0.1. 
These four observations alone place strong constraints on 
the physical properties of L1689B, as we show in the next 
section. 

In particular we wish to explore two different possibili- 
ties. When low contrast is observed in molecular maps, and 
the brightness temperature of different transitions is similar 
- as is true in this case - it is often assumed that this implies 
that the lines are optically thick. However another possibil- 
ity is that the optical depth of each line is constant across 
the map, a situation which can arise if: (i) temperature gra- 
dients affect the fraction of molecules in each level so as to 
maintain a relatively constant column density of molecules 



in upper excited levels across the core; or (ii) the molecu- 
lar abundance changes so as to keep the column averaged 
abundance of CO itself constant. The ratio of the line inten- 
sities may also appear equal in these situations - if T^x ~9K 
the optical depth of the (1— >-0) line is equal to the optical 
depth of the (2— >1) line, and if Tex ~23K the optical depth 
of the (2— >1) line is similar to that of the (3^2) line. For 
the remainder of this paper we use a parameterised model 
of the core, in combination with a radiative transfer code to 
investigate these effects and others in detail, and we argue 
that a molecular abundance drop is the the most plausible 
explanation for the appearance of L1689B. This abundance 
drop is most likely to be due to the freeze-out of molecules 
onto the dust grains, giving further evidence of the prestellar 
nature of L1689B. 



3 MODEL DESCRIPTION 
3.1 Model geometry of L1689B 

We have chosen a model which allows certain key physi- 
cal characteristics to be represented in terms of simple an- 
alytical functions, and which at the same time incorporates 
qualitatively various physical models suggested, whilst lim- 
iting the number of free parameters. In this way we can pa- 
rameterise the observed appearance of L1689B and discover 
how the observations constrain these parameters. We can 
then hope to discover which (if any) of the theoretical mod- 
els of protostar formation best describes the observations of 
L1689B. 

The model parameters are shown schematically in Fig- 
ure 4. We first make the simplification that L1689B is spher- 
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Figure 4. The model core for L1689B, partly based on the find- 
ings of Paper II. The core has an outer envelope with a density 
profile which falls as r~^. The outer envelope is isothermal and 
is assumed to have the CO abundance typical of the Ophiuchus 
cloud. The inner core region has a more flattened density profile 
as suggested by the mm/submm continuum observations. The 
inner core is not assumed to be isothermal, and the CO abun- 
dance can also vary. In order to keep all parameters finite, the 
core has an inner radius of 0.004pc, well below any current single 
dish telescope resolution, within which all physical conditions are 
homogeneous. 

ically symmetric, which significantly simplifies the radiative 
transfer analysis. Given that the ellipticity of L1689B ap- 
pears to be low (Paper II), we feel this justifies the consider- 
able saving of coding and cpu time. The model core consists 
of an isothermal envelope of outer radius 0.07pc surrounding 
an inner core of radius 0.02pc (c.f. Paper II). The tempera- 
ture and CO abundance are kept constant in the envelope, 
and the density varies as r~^. Within the inner core, all three 
physical parameters of density, temperature and abundance 
are allowed to vary according to a power-law dependence, 
as shown in Figure 4. The radius at which a break in the 
power law density profile is observed (the break between the 
inner core and outer envelope) will also be referred to as the 
critical radius (c.f. R/jat in Paper II). To ensure that the 
parameter values remain finite we also set an inner radius 
of the core equal to 0.004pc (well within the resolution of 
JCMT or IRAM), which does not affect the results. 

The inner core density profile, temperature profile and 
CO abundance profile power-laws are given by , r°"^, 
and r"^ respectively. Hence there are 5 free parameters in 
the model: the central density; the outer temperature; and 
the value of the three power-law indices. The critical radius 
between core and envelope is set by the mm continuum data 
(Paper II). The outer abundance is set to the canonical C^®0 
abundance for Ophiuchus. 

3.2 The A-Iteration Code 

We used a previously tried and tested radiative trans- 
fer code, known as the Stenholm Code (Stenholm 1977, 



Matthews 1986), which uses a A- iteration method to solve 
for the population levels and produce model spectra output. 
The code models a spherically symmetric molecular cloud 
as a set of shells, each with uniform physical conditions. 
For each shell it requires that the H2 density, the H2 tem- 
perature and the abundance of C^^O (with respect to H2) 
be specified. Because molecular line widths are generally 
greater than predicted from simple thermal line broaden- 
ing, a micro-turbulent component of velocity, 5v, can also 
be specified for each shell. 

The code starts by setting the population levels of the 
C^^O rotationally excited states to be consistent with the 
Boltzmann distribution. The equation of radiative transfer 
is then numerically solved to obtain 1^, the radiation field 
in each shell. Using this radiation field the rate equation is 
then solved to give a new, modified set of population levels. 
The procedure of calculating nj (the population levels of the 
C^^O) followed by I^ (the radiation field) is repeated until 
a stable solution is approached and convergence is achieved. 
The code can then solve the radiative transfer equation for 
lines of sight through the cloud and simulate observed line 
profiles by smoothing with a Gaussian spatial profile similar 
to the beam profile of the JCMT. The H2 collision rates for 
the model were obtained from Flower & Launey (1985), who 
calculated the collision rates for para- 1/2 up to the J = 11 
rotational level, for temperatures between 10 and 250K, and 
for ortho-/i'2 up to the J = 6 level between 10 and lOOK. 

Throughout the simulations (unless stated otherwise) 
the model core representing L1689B had the following prop- 
erties: The model core consisted of 20 logarithmically spaced 
shells (giving increased resolution towards the centre of the 
core), and a critical radius rcrit = 0.02pc marking the bound- 
ary between the core and envelope. The micro-turbulent ve- 
locity component for each shell throughout the cloud was 
given by 5v(r) oc r" '', up to a maximum 0.8 kms^^ at the 
outside. This was chosen to be consistent with the observa- 
tions summarised by Larson (1981) and to fit the observed 
line widths. The ortho- to para- 7^2 ratio was fixed at 1. 

The following section investigates the effect of varying 
the free parameters in the model of L1689B: qt; ftp; a^; 
and Toutor- The central density, pc, was normalized in the 
simulations so as to reproduce the central brightness of the 
core in the (2-^1) transition. 

4 MODEL RESULTS 
4.1 Isothermal Models 

The simplest models of cloud cores assume a constant tem- 
perature. So we first investigate how the predicted appear- 
ance of L1689B is dependent on the inner core radial density 
profile p oc r~"'' if a single temperature, T {=T outer), exists 
throughout the core. In other words, we set both ar and 
a-^ to zero. Hence, the predicted appearance for the set of 
cores with temperatures between 10 and 25K and Op be- 
tween and 2 was calculated. By the nature of this form of 
modelling, we can run whole families of models and produce 
large output data-sets. We then illustrate these models as a 
series of contour plots in the two-dimensional phase space 
of the two parameters we are varying. 

Figure 5 shows two contour plots that illustrate the ef- 
fects of varying Qp and T. Figure 5(a) shows contours of the 
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T(32 arcsecs)/T(centre) 




Figure 5. Two contour plots illustrating the effects of varying Qp and T, whilst setting qx = ctx ~ (i-^- isothermal models), (a) 
Contours of the values of the ratio of Tjnb(J=3— >2) to Tiiib(J=2^1). The area of this plot that is consistent with the data shown in 
Figures 1 & 2 is shaded, (b) Contours of the values of the ratio of Tj^b at a radius of 32 arcsec from the centre to Tj^i, at core centre for 
the C'^^O (J=2— >1) transition. The value observed in the L1689B data for this ratio was l.OiO.l. No physically realistic value of Up can 
reproduce this value of the ratio. 



values of the ratio of Tn,b(J=3^2) to Tn,b(J=2^1). The 
area of this plot that is consistent with the data shown in 
Figures 1 & 2 is shaded. Figure 5(b) shows contours of the 
values of the ratio of T^b at a radius of 32 arcsec from the 
centre to T^b at core centre for the C^**0 (J=2-^l) tran- 
sition. As noted in section 2 above, the value observed in 
the data of L1689B for this ratio was l.OiO.l. It can be 
seen from Figure 5(b) that for T>13K this value increases 
as ttp decreases and that the brightness profile in this area is 
mainly dependent on the density profile of the core. However 
for T<13K the brightness profile becomes strongly depen- 
dent on T. This is because the C^^O transitions become 
increasingly optically thick. 

The reason why the (J=3^2) to (J=2— »1) brightness 
ratio of the model core increases with both T and Qp is be- 
cause the upper rotational levels of C^*0 become more pop- 
ulated as the temperature and density increase. It should be 
noted that the singular isothermal sphere model represents 
the upper border of Figure 5(b). This is the part of the plot 
furthest from consistency with the data for T>13K. This 
agrees with the conclusions of Papers I & II that L1689B 
does not represent the singular isothermal sphere initial con- 
ditions for star formation. 

Although at first sight setting T=10K seems to offer a 
possible explanation for L1689B's brightness profile, another 
effect led us to reject this as a possible model for L1689B. 
At lOK the model has very high optical depths and the 
predicted line profiles are significantly wider and show an 
extended 'flat top'. This is illustrated in Figure 6. The cen- 
tral C^*0 line profile is shown as observed (the data points) 
and as predicted by the model ap=0, T=12K (solid line) 
and ap=0, T=10K (dashed line). It is clear that the T=10K 
model does not agree with the data. We therefore need to 
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Figure 6. Comparison of observed C^^O (J=2^1) line (data 
points) with two isothermal models. The first model (solid line) 
has T=12K, the second (dashed) has T=10K. The second clearly 
shows a flat top due to high optical depths, and is not consistent 
with the data. 



find some other mechanism to explain the low contrast ob- 
served in the data. 



4.2 Varying temperature models 

We next investigated a set of models in which the outer en- 
velope was at a temperature of 20K, but in the inner region 
the temperature profile falls as T oc r""^ to a minimum tem- 
perature of lOK. This means that, when ot is greater than 
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Figure 7. Two contour plots illustrating the effects of varying ap and aj- simultaneously, (a) Contours of the values of the ratio of 
Tmb(J=3— >2) to Tnib(J=2— >!)• The shaded area is consistent with the L1689B data, (b) Contours of the values of the ratio of Tnib a-t a 
radius of 32 arcsec from the centre to Tj^b at core centre for the C^*0 (J=2— >1) transition. The shaded area is consistent with the data 
from Figures 1 & 2. Note that the two shaded areas overlap slightly. 



approximately 0.5, an inner lOK region is created, whose 
size increases as ot rises. When ot is 2.0, the radius of this 
region is 0.015pc. 

Figure 7 shows two contour plots that illustrate the ef- 
fects of varying Qp and qt simultaneously. Figure 7(a) shows 
contours of the values of the ratio of T^b ( J=3^2) to T^b 
(J=2-^l). The area of this plot that is consistent with the 
L1689B data is shaded. 

Figure 7(b) shows contours of the values of the ratio of 
Tmb at a radius of 32 arcsec from the centre to Tmb at core 
centre for the C^*0 (J=2— >1) transition. As noted in section 
2 above, the value observed in the data of L1689B for this 
ratio was l.OiO.l. Once again the part of this plot that is 
consistent with the data from Figures 1 & 2 is shown as a 
shaded area on Figure 7(b). The fact that figure 7(b) shows a 
much larger variation than Figure 7(a) illustrates the impor- 
tance of taking observations at more than one radius if one 
wishes to truly characterize cloud cores, and shows clearly 
that single spectra at the centre of a cloud do not constrain 
its parameters. 

This time there is a set of values of Up and ar that 
is consistent with all of the C^^O data. This approximately 
corresponds to 1.1 < ar < 2 and < Qp < 0.4. However, at 
the la level, we can also say that if Qp < 0.1, then qt must 
by < 1.5. These results reveal that a decrease in temperature 
towards the centre of the core does indeed help to explain 
the C^'^O appearance of L1689B. In fact a model which has 
a flat inner density profile and a very steep temperature 
drop from 20K at a radius of 0.02pc, to lOK at a radius of 
0.015pc, is consistent with our observations. 




Figure 8. Sketch of the assumed observing geometry of L1689B. 



4.3 Continuum data revisited 

A 1.3-mm continuum map of L1689B was published in Pa- 
per II. In this section we incorporate those observations to 
see what further constraints we can place on our model of 
L1689B. Remembering that the millimetre continuum arises 
from thermal emission from dust grains, the flux density in 
this situation (e.g. Hildebrand 1983) is given by: 

= B4T)n(l-e-"-), (2) 
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T(32 arcsecs)/T(c = ntr=) 




0.5 1 1.5 2 



Figure 9. Contours of the values of the ratio of Tj^i, at a radius 
of 32 arcsec from the centre to T^^, at core centre for the C^*0 
(J=2— >1) transition. The shaded area at lower right is consistent 
with the data from Figures 1 & 2. The shaded strip towards the 
upper left indicates the constraints placed on ax and Op by the 
continuum data. Note that the two shaded areas do not overlap. 



If the density profile p(r) and the temperature profile 
T(r) are known, where r is the radial distance from the cloud 
centre, then the observed surface brightness I(R) is given 
by the equation of radiative transfer in the assumption of 
optically thin emission: 

p{r)B{T[r])dx. (4) 

-XQ 

Making the substitutions R = r cos e and dx = [Rde/cos^e] 
we derive: 

I{R)(xR / p{R/ cos 6) B{T[R/ cose]) —.(5) 

For the model of L1689B, the density profile in the outer 
envelope is represented by the power law p oc r~^. Hence, 
for R > rcrit (i.e. for observations of the outer envelope) we 
have: 

I{R) oc R~^ / (cos6')"''d6i, (6) 

J-7r/2 

where the limits in the integration have been set using the 
approximation that the core extends to infinity. This integral 
is independent of R and therefore we predict the I{R) oc 
R~^ behaviour in the envelope of L1689B that was actually 
observed in the continuum in Paper II. 

However, for the inner core the equation for the pro- 
jected surface brightness becomes more complicated: 



where Bi,{T) is the black-body function at temperature T, 
Q. is the source solid angle and is the optical depth. 
This form is often referred to as a grey-body function. The 
black-body function can often be simplified to the so-called 
Rayleigh- Jeans approximation in cases where {hv/kT) « 
1. However, at 1.3 mm, {hv/k) = ICSif, which is within a 
factor of 2 of our fitted temperature for the the L1689B core 
(see above) . Hence, in this case the Rayleigh- Jeans approxi- 
mation may introduce errors, so we here outline a somewhat 
more rigorous approach. 

Figure 8 shows our assumed observing geometry. We 
still make the assumption that the dust properties, such as 
grain size distribution, abundance with respect to H2 and 
emissivity, are constant throughout the cloud. Then we see 
that the intensity arising from an element of the cloud of 
length dx at a given frequency is: 

dl{x) oc p(x)B{T{x))dx, (3) 

where p{x) is the density profile along the line of sight and 
T{x) is the temperature profile of the dust along the line of 
sight - see Figure 8. 

The simple form often used for power-law radial profiles, 
is that if T{r) oc r~', p{r) oc r~^, and the fiux density 
Sv{9) oc then the indices are related by the simple 

expression m = p + q— 1 (Casali 1986, Adams 1991). This is 
derived by assuming the depth of the cloud remains constant 
across the entire core (e.g. a slab or infinite cloud). However, 
this form only holds for the Rayleigh- Jeans approximation, 
and for the case where the radial density and temperature 
each follows a single power-law, and there is no break in the 
power laws, such as we have in our model - see Figure 4. 
We derive below a more rigorous set of relations. 



I^R)^^(2^(lcos~\.))+ (7) 

^ In 23 "^1 

JO cos-^Q ' 

where we have used the substitution — R/ rcrit, and Tcnv 
is the envelope temperature. The first term of this equation 
is the contribution to the signal from the envelope, while the 
second term is the contribution from the core. This equation 
can be solved numerically for any value of Op or ar, to 
predict the brightness profile I(R). 

Numerical solution of equation 7 in the isothermal case 
allows us to compare the continuum radial fiux density varia- 
tion observed in Paper II, with that predicted by the param- 
eterised model we have presented here. We find that models 
with p oc r-°-^^°-'^ (i.e. Op = 0.5 ± 0.1) best match the 
data presented in Paper II. A more general, non-isothermal, 
solution to the continuum appearance of L1689B can be ex- 
pressed approximately as: 

Qp - 1.5QT = 0.5±0.1. (8) 

This solution was found by numerically solving equation 7 
for all values of Op between and 2, and ot between and 
2, and then finding which range of values give solutions best 
matching the data in Paper II. This is a new and different 
constraint to that discovered in the above modelling of the 
C^^O data. 

Figure 9 illustrates this additional constraint introduced 
by the continuum data by reproducing the contours from 
Figure 9(b) - the values of the ratio of Tmb at a radius 
of 32 arcsec from the centre, to Tmb at core centre for the 
C^^O (J=2— >1) transition - as a function of varying ap and 
Qt- Once again the part of this plot that is consistent with 
the CO data from Figures 1 & 2 is shown as a shaded area 
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Figure 10. Two contour plots illustrating the effects of varying Op and simultaneously, (a) Contours of the values of the ratio of 
Tinb(J=3— >2) to T^b(J=2— >1). The shaded area is consistent with the L1689B CO data, (b) Contours of the values of the ratio of Tnib 
at a radius of 32 arcsec from the centre to Ty^h core centre for the C^*0 {J=2— >1) transition. The shaded area at the lower right 
hand corner is consistent with the data from Figures 1 & 2. The horizontal shaded strip indicates the constraints placed on ap by the 
continuum data. 



on Figure 7(b). However, now we have added an additional 
shaded strip at the upper left of the plot illustrating the 
constraint placed by the continuum data, as described by 
equation 8. We see that the two shaded areas do not overlap, 
indicating that there is no simultaneous fit to the CO data 
and the continuum data. Hence, a more complex model is 
required to simultaneously fit all of the data. 

4.4 Varying the CO abundance 

The simplest way to try to reconcile the CO data to the 
continuum data is to vary the abundance of gas-phase CO 
relative to H2 . This is predicted to occur in regions of high 
density, as the CO freezes out onto the dust grains, and 
it has been observed qualitatively in many regions on many 
occasions (e.g. Mezger et al. 1990; Gibb & Little 1998; Ward- 
Thompson et al. 2000). It is thought to occur at tempera- 
tures less than 17K (Nakagawa 1980, Bergin & Langer 1997). 
So we ran a set of models to investigate the efi'ect of decreas- 
ing the CO abundance towards the centre of the core. These 
models had a single temperature of 20K, an inner density 
gradient as before, of p oc r~°'p , and a CO abundance frac- 
tion X oc r"^ . We then investigated the regime < cip < 2 
and < Qx < 2. 

Figure 10 shows two contour plots that illustrate the 
effects of varying Up and simultaneously. Figure 10(a) 
shows contours of the values of the ratio of Tmb ( J=3^2) to 
Tmb(J=2^1). The area of this plot that is consistent with 
the CO data is shaded. 

Figure 10(b) shows contours of the values of the ratio 
of Tmb at a radius of 32 arcsec from the centre to Tmb at 
core centre for the C^*0 (J=2— >1) transition. As noted in 
section 2 above, the value observed in the data of L1689B 
for this ratio was l.OiO.l. Once again the part of this plot 



that is consistent with the data from Figures 1 & 2 is shown 
as a shaded area in the lower right hand corner of Figure 
10(b). 

Including the continuum data in this case, with qt = 0, 
equation 10 simplifies to ap — 0.5 ± 0.1. Figure 10(b) il- 
lustrates this additional constraint as an horizontal shaded 
strip in the lower half of the plot. Again we see that the 
two shaded areas do not overlap, indicating that there is no 
simultaneous formal fit to the CO data and the continuum 
data. However, this is only at the la level, and we find that 
the shaded regions would overlap at the 2a level for the val- 
ues ~ 1.9 — 2 and ap ^ 0.4. This abundance effect is quite 
extreme - it implies a 95% drop in C'^^O abundance from 
edge to centre - although not unprecedented. Gibb and Lit- 
tle (1998) found a similar reduction in the C^**0 abundance 
in HH23-26. 

Furthermore, if there is such a large amount of freeze- 
out taking place, this will increase the average grain size to- 
wards the core centre and consequently the dust grain mass 
opacity. In this case we would over-estimate the central den- 
sity profile exponent, Up derived from the continuum data. 
This would then act to shift the horizontal shaded region in 
Figure 10(b) downwards, causing a genuine overlap region 
at the 1-a level. Thus we have a fit to the data with 2 
and ftp < 0.4. The spectra produced by the model were 
compared with those observed and there was good agree- 
ment within the errors. This is an isothermal fit. However 
we cannot rule out a slight temperature gradient, depend- 
ing upon how much the dust opacity parameters that are 
not well constrained can vary. 

We ran this set of simulations for several temperatures 
from 10 to 25K. We found that the predicted surface bright- 
ness did not vary significantly with T, even at very low tem- 
peratures. However, the predicted ratio of C^*0 (J=3-^2) to 
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C^^O (J=2-^l) did fall with T, and for models with T<14K 
the predicted values were too low to match our observations. 



4.5 Other factors 

There are other factors than those discussed above that may 
affect the conclusions of the modelling. One such factor is 
that the ortho- to para-H2 ratio in these cores is not well 
defined. Since ortho-H2 has a higher cross section for colli- 
sions with CO, the ratio of ortho- to para-H2 may affect the 
conclusions drawn. It was found that by increasing the frac- 
tion of ortho-H[2, the C^^O centre-to-edge ratio increased. 
However, it hardly affected the ratio of the the (J=3— >2) 
to (J=2— >1) transitions. Varying the ortho- to para- ratio 
altered the centre-to-edge ratio by 0.1 at most. 

Similarly, the precise form of the microturbulence pro- 
file may infiuence the appearance of the core. For exam- 
ple, it was found that when using a constant 5v microtur- 
bulence profile, the brightness gradient decreased slightly - 
T(32 arcsec)/T(centre) increased by ~ 0.1. However, neither 
of these two factors affected the appearance of the L1689B 
data strongly enough to affect the conclusions drawn above. 

The geometry of the core as a whole may also change 
the exact form of the solutions. In particular, any departure 
from spherical symmetry would have an effect. Likewise, any 
variations from homogeneity - i.e. dumpiness within the 
core - is also likely to affect the predicted appearance of the 
core in mm/submm observations. It is difficult to quantify 
these effects, and further data would be needed to constrain 
these extra parameters. 

We also investigated the possible effect of having a small 
temperature gradient in the envelope, or allowing the density 
gradient in this region to fall less steeply (p oc n"^'^). This 
was found to only vary the ratio of T(32 arcsec)/T (centre) 
by ~ 0.05. Checks to ensure that the exact values of the 
inner radius, and the cloud size did not significantly affect 
our results were also made. 



5 CONCLUSIONS 

In this paper we have presented new C^**0 (J=3-*2) and 
(J=2— >1) data of the pre-stellar core L1689B. We have used 
a spherically symmetric radiative transfer code to model 
these data in terms of three parameters - the gradients in 
temperature, density and C^^O abundance. In the three di- 
mensional parameter space defined by the three exponents, 
Qp, qt, and a^, we have found solutions consistent with the 
data. There is a region defined by Up — 1.5qt = 0.5 ± 0.1, 
which is consistent with the mm continuum observations of 
L1689B. 

There is also a set of solutions that can simultaneously 
predict the molecular CO appearance and the mm contin- 
uum data, implying ap < 0.4, ar = (although this is 
not well constrained) and a,^ = 2. Hence, we have shown 
that there is freeze-out of CO towards the centre, and we 
have constrained the radial density and temperature profiles. 
Thus we have shown that the combination of mm/submm 
spectral line and continuum data, with a rigorous treatment 
of radiative transfer, is a powerful method of investigating 
the physics of pre-stellar evolution. 
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